import numpy as np

#from scipy import optimize
#import scientific_utils as su
#from matplotlib.pyplot import *
#import spectrometer as sp
#reload(su)
#import copy

def propagate(source,target):
	E = source.wave.E
	dE = source.wave.dE
	k = source.wave.k
	Q = source.curve.point
	dQ = source.curve.width
	n = source.curve.normal
	P = target.curve.point
	m = target.curve.normal
	
	P.shape += (1,)
	Q.shape += (1,)
	P.transpose(0,2,1)
	R = P-Q
	r = np.sqrt(np.sum(R**2,0))
	cosrn = sum(R*n,0)/r
	cosrm = sum(r*m,0)/r
		
	U = np.sum(np.exp(1j*k*r)*np.sqrt(k/r)*(1j*E*cosrn/2-dE)*dQ,0)
	dU = 1j*k/2*U*cosrm

	target.wave.k = k
	target.wave.E = U
	target.wave.dE = dU

#def ApperturesPropagate(SourceWave,k,SourceCurve,TargetCurve):
#	SourceCurveType = type(SourceCurve)
#	if SourceCurveType == 'tuple':
#		for i in SourceCurve:
#			propagate(SourceWave,SourceCurve,TargetCurve)
#	elif SourceCurveType == 'LineSegment':
#		pass
#	elif SourceCurveType == 'Curve':
#	
#	TargetCurveType = type(TargetCurve)	
#	if TargetCurveType == 'tuple':
#		for i in TargetCurve:
#			propagate(SourceWave, SourceCurve, TargetCurve)
#	
#	E = SourceWave[0]
#	dE = SourceWave[1]
#	Q = SourceCurve[0:1]
#	n = SourceCurve.normals()
#	dQ = SourceCurve[3]
#	TargetCurve[0:1]
#	for p in TargetCurve:
#		P = p[0:1]
#		m = p.normals()
#		r = P-Q
#		cosrn = inner(r,n)/length(r)
#		cosrm = inner(r,m)/length(r)
#		r = length(r)
#		U = sum(exp(1j*k*r)*sqrt(k/r)*(1j*E*cosrn/2-dE)*dQ)
#	dU = 1j*k/2*U*cosrm
#	return vstack((U,dU))			
#
#def ApperturesPropagate		
#
#def propagate(SourceWave,LineSegment,TargetCurve):
#	
#class WavePropagator:
#	def __init__(self,spectrometer):
#		self.E=0
#		self.keff=0
#		self.spectrometer=spectrometer
#		self.spectrometer.grating.E = [0]
#		self.setAperturePoints(spectrometer.input,100)
#		self.setCosMode(spectrometer.input)
#		self.spectrometer.input.E = self.spectrometer.input.mode
#		self.eval_r1=True
#		self.eval_r2=True
#
#	def setCosMode(self,aperture):
#		aperture.mode = cos(pi*aperture.x/aperture.width)*sqrt(2/aperture.width)
#		
#	def setAperturePoints(self,aperture,n):
#		aperture.n = n
#		self.spectrometer.setCenter(aperture)
#		aperture.eval_r2 = True
#		angle = aperture.angle + pi/2
#		aperture.x = linspace(-aperture.width/2,aperture.width/2,n)
#		aperture.points = transpose(vstack((cos(angle)*aperture.x,sin(angle)*aperture.x)))+aperture.center		
#		
#	def propagateToGrating(self,keff):
#		if self.keff==keff:
#			return
#		if self.eval_r1:
#			r1 = self.spectrometer.grating.points - self.spectrometer.input.center
#			self.alpha = arctan2(r1[0],r1[1])-self.spectrometer.input.angle
#			self.alphaI = (su.angle(-r1)-self.spectrometer.grating.angles)
#			r1 = sqrt(sum(r1**2,1))
#		if len(self.spectrometer.grating.E) != self.spectrometer.grating.n:
#			self.spectrometer.grating.E=zeros(self.spectrometer.grating.n,dtype='complex')
#		for i in range(self.spectrometer.grating.n):
#			dr1 = sin(self.alpha[i])*self.spectrometer.input.x
#			self.spectrometer.grating.E[i] = sum(self.spectrometer.input.mode*exp(-1j*keff*dr1))
#		dx = self.spectrometer.input.x[1]-self.spectrometer.input.x[0]
#		self.spectrometer.grating.E = self.spectrometer.grating.E*(1+cos(self.alpha))*exp(-1j*keff*r1)*dx
#		self.spectrometer.grating.E = self.spectrometer.grating.E/(2*sqrt(2*pi*r1/keff))
#		self.keff = keff
#			
#	def propagateTo(self,target):
#		if not hasattr(target,'E'):
#			target.E = zeros(target.n,dtype='complex')
#		if len(target.E) != target.n:
#			target.E = zeros(target.n,dtype='complex')
#		if target.eval_r2:
#			target.eval_r2 = False
#			target.r2 = zeros((target.n,self.spectrometer.grating.n))
#			target.alphaD = zeros((target.n,self.spectrometer.grating.n))
#			target.sin_diff = zeros((target.n,self.spectrometer.grating.n))
#			target.cos_sum  = zeros((target.n,self.spectrometer.grating.n))
#			for i in range(target.n):
#				r  = target.points[i]-self.spectrometer.grating.points
#				target.r2[i] = sqrt(sum(r**2,1))
#				target.alphaD[i] = (self.spectrometer.grating.angles-su.angle(r))
#				target.sin_diff[i] = (sin(self.alphaI)-sin(target.alphaD[i]))/2
#				target.cos_sum[i] = (cos(self.alphaI)+cos(target.alphaD[i]))
#		for i in range(target.n):
#			EfromEachFacet = sinc(self.keff*self.spectrometer.grating.widths*target.sin_diff[i]/pi)
#			EfromEachFacet = EfromEachFacet*self.spectrometer.grating.widths
#			EfromEachFacet = EfromEachFacet*exp(-1j*self.keff*target.r2[i])/sqrt(target.r2[i])
#			EfromEachFacet = EfromEachFacet*self.spectrometer.grating.E*target.cos_sum[i]
#			target.E[i] = sum(EfromEachFacet)*sqrt(self.keff/8/pi)
#			
#	def fractionCoupledInto(self,waveguide):
#		## Calculate the internal product using the simpson rule.
#		if waveguide.n%2==0:
#			print 'Number of points need to be odd.'
#			return -1
#		h = waveguide.x[1]-waveguide.x[0]
#		weights = [1, 4] + (waveguide.n-3)/2*[2, 4] + [1]
#		return (sum(weights*waveguide.mode*abs(waveguide.E))*h/3)**2
#		
#	def transmissionSpectrum(self,keffs,waveguide):
#		spectrum = zeros(len(keffs))
#		for i in range(len(keffs)):
#			self.propagateToGrating(keffs[i])
#			self.propagateTo(waveguide)
#			spectrum[i] = self.fractionCoupledInto(waveguide)
#		return spectrum
#	
#	def findPeakTransmissioAngleAt(self,keff,aperture):
#		self.propagateToGrating(keff)
#		N = su.fwhm(abs(self.spectrometer.grating.E))
#		deltaLambdaEff = 2*pi/(self.spectrometer.order*N*keff)
#		a = self.spectrometer.grating.pitch
#		m = self.spectrometer.order
#		lambdaEff = 2*pi/keff
#		inputAngle = self.spectrometer.input.angle
#		startAngle = arcsin(m/a*(lambdaEff-deltaLambdaEff/2)+sin(inputAngle))
#		endAngle   = arcsin(m/a*(lambdaEff+deltaLambdaEff/2)+sin(inputAngle))
#		def f(angle):
#			self.spectrometer.changeApertureAngle(aperture,angle)
#			self.propagateTo(aperture)
#
#			return -self.fractionCoupledInto(aperture)
#		return optimize.golden(f,brack = (startAngle,endAngle), tol=1e-6,full_output=True)
#		
#	def optimizeWidth(self,aperture):
#		self.propagateToGrating(aperture.keff)
#		def f(width):
#			aperture.width = width
#			self.setAperturePoints(aperture,21)
#			self.setCosMode(aperture)
#			self.propagateTo(aperture)
#			a = self.fractionCoupledInto(aperture)
#			print 'fraction =', a
#			return -self.fractionCoupledInto(aperture)
#		return optimize.golden(f,brack = (100e-9,10e-6), tol=1e-6)		
#
#	def setOutputArray(self,keffs):
#		aperture = sp.Aperture()
#		aperture.width = 1.4e-6
#		aperture.angle = 30*pi/180
#		self.spectrometer.setCenter(aperture)
#		self.setAperturePoints(aperture,51)
#		self.setCosMode(aperture)		
#		self.spectrometer.outputs = []
#		i=0
#		for keff in keffs:
#			ape = copy.deepcopy(aperture)
#			ape.keff = keff
#			self.findPeakTransmissioAngleAt(keff,ape)
#			
#			#self.optimizeWidth(ape)
#			self.spectrometer.outputs.append(ape)
#			i+=1
#			print i,' done' 
#		
#		
#				
#	def plotInputField(self):
#		pass
#	
#	def plotGratingField(self):
#		figure()
#		plot(abs(self.spectrometer.grating.E))
